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ABSTRACT 

Based on numerical simulations of galaxy mergers that incorporate black hole growth, we predict the faint 
end slope of the quasar luminosity function (QLF) and its evolution with redshift. Our simulations have yielded 
a new model for quasar lifetimes and light curves where the lifetime depends on both the instantaneous and peak 
luminosities of the quasar. This description motivates a new interpretation of the quasar luminosity function 
in which the bright end of the QLF consists of quasars radiating near their peak luminosities, but the faint 
end is mostly made up of brighter peak luminosity quasars seen in less luminous phases of evolution. The 
faint-end QLF slope is then determined by the faint-end slope of the quasar lifetime for quasars with peak 
luminosities near the observed break. We determine this slope from the quasar lifetime as a function of peak 
luminosity, based on a large set of simulations which encompass a wide variety of host galaxy, merger, black 
hole, and interstellar gas properties. Brighter peak luminosity (higher black hole mass) systems undergo more 
violent evolution, and expel and heat gas more rapidly in the final stages of quasar evolution, resulting in a 
flatter faint-end slope (as these objects fall below the observed break in the QLF more rapidly). Therefore, 
as the QLF break luminosity moves to higher luminosities with increasing redshift, implying a larger typical 
quasar peak luminosity, the faint-end QLF slope flattens. From our determined quasar lifetimes as a function 
of peak luminosity and this interpretation of the QLF, we predict the faint-end QLF slope and its evolution with 
redshift in good agreement with observations. Although black holes grow in an anti-hierarchical manner (with 
lower-mass black holes formed primarily at lower redshifts), the observed change in slope and differential or 
"luminosity dependent density evolution" in the QLF is completely determined by the non-trivial, luminosity- 
dependent quasar lifetime and physics of quasar feedback, and not by changes in the shape of the underlying 
peak luminosity or active black hole mass distributions. 

Subject headings: quasars: general — galaxies: active — galaxies: evolution — cosmology: theory 



1. INTRODUCTION 

The shape and evolution of the quasar luminosity func- 
tion (QLF) is fundamental to cosmology, constraining theo- 
ries of galaxy and supermassive black hole (BH) formation, 
accretion models, the X-ray, UV, and infrared backgrounds, 
the role of galaxy mergers and interactions, and reioniza- 
tion. Until recently, it has not been possible to reliably mea- 
sure the faint-end slope of the QLF even at low redshifts, 
but this has begun to change with the advent of large, uni- 
formly selected quasar samples in surveys such as the SDSS 
and 2dF which now increasingly probe the QLF to fainter lu- 
minosities. Furthermore, a growing number of observations 
at different redshifts, in radio, optical, and soft and hard X- 
ray bands, have suggested that the faint end slo pe evolves , 
becoming flatter at higher redshift (e .g. Page et alJ 119971 
Mivaii et al."2000', 2001; "La Franca et al. 2002; Cowie et al.' 
2003; Ueda et al. 2003; Fiore et al. 2003; Huntetal. 2004; 
Cirasuolo e t al. 2005; Hasinger et al. 2005). 

This evolution, parameterized either as direct evolution in 
the faint-end slope or as "luminosity-dependent density evo- 
lution" (LDDE), has been the subject of much speculation, as 
it implies that the density of lower-luminosity quasars peaks 
at lower redshift. In traditional models of the quasar life- 
time and light curve, this evolution must be directly related 
to the distribution of quasar hosts, implying e.g. significant 
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and rapid evolution in the shape of the distribution of host 
galaxy masses, which cannot be accounted for in either semi- 
analytical models or numerical simulations and is not con- 
sistent with a wide range of galaxy observations. Although 
models which adopt these simplified prescriptions for quasar 
evolution have had success predicting the evolution of the 
bright end of the QL F (e.g. iKauffmann & Haehnelt 200(i 
'Wvithe & Loeb'2003i, they do not predict the faint-end slope 
or its evolution, and as such cannot be reliably extrapolated 
to low luminosities or to redshifts where the slope is unde- 
termined. Observations of the high-redshift, faint end quasar 
luminosity function slope are highly uncertain, and no large, 
uniformly selected samples yet exist which measure the faint- 
end slope at both low (z < 1) and high (z > 3) redshifts, and 
therefore theoretical models of the faint-end slope are espe- 
cially important. This high-z, faint-end slope is a critical 
quantity in determining the early formation history of black 
holes and, especially, their contribution to reionization, as 
well as possible connections between quasars and e.g. the 
low-luminosity Seyferts seen at z ~ 0. 

Without a more sophisticated model of quasar evolution, 
models which have attempted to reconcile observations of 
evolution in the faint-end QLF slope and BH populations 
(e.g., Merloni 2004) have been forced to fit to the entire QLF 
and BH mass distribution to essentially arbitrary distributions 
of lifetimes/duty cycles and accretion rates as a function of 
redshift. Still, this phenomenological modeling has eluci- 
dated the anti-hierarchical nature of BH growth, with smaller- 
mass BHs formed at lower redshift as an implication of this 
evolution in the QLF. But an actual prediction of the faint-end 
slope requires a more detailed model for both quasar lifetimes 
and the QLF 
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Recently, BH growth and feedback have been incorporated 
into numerical simulations of galaxy mergers (Springel et al. 
2005a,b, Springel 2005, Di Matteo et al. 2005). In these sim- 
ulations, gravitational torques drive inflows of gas into the 
nuclei of merging galaxies (e.g. Barnes & Hernquist 1991, 
1996), triggering starbursts (e.g. Mihos & Hernquist 1996) 
and feeding the growth of central supermassive BHs (Di Mat- 
teo et al. 2005). As the BHs accrete, some of the radiated 
energy couples to the surrounding gas, and the growth even- 
tually stalls when this feedback energy is sufficient to un- 
bind the reservoir of gas. These simulations elucidate the 
connection between galaxy evolution, the formation of super- 
massive BHs, and the self-regulated nature of quasar activity, 
and provide quantitative predictions which agree well with 
observations of, e.g., the Mbh - c relation (Di Matteo et al. 
2005, Robertson et al. 2005), quasar lifetimes (Hopkins et al. 
2005a,b), and the QLF in various wavebands (Hopkins et al. 
2005c,d,e). These simulations provide a fully self-consistent, 
quantitative prediction of the light curve of quasar activity and 
the resulting quasar "lifetime" - i.e., the amount of time that 
a quasar spends at a given luminosity - and its dependence on 
the properties of the merging galaxies. 

The resulting quasar light curves imply a qualitatively 
different picture of t he QLF than previously considered 
jHopkins et al.' "2005c'). Specifically, in the model of Hop- 
kins et al. (2005a-f), quasars evolve rapidly and their life- 
time depends on both their instantaneous and peak luminosi- 
ties. Critically, the quasar lifetime in this model is longer at 
lower luminosities; i.e. a given quasar spends more time (and 
is more likely to be observed) at a luminosity well below its 
peak luminosity. This differs from previous models which 
generally assume that quasars radiate at a fixed luminosity for 
some characteristic lifetime, or adopt simplified exponential 
light curves. In our picture, the bright end of the QLF consists 
of quasars growing rapidly and radiating near their peak lumi- 
nosities, while the faint end consists of quasars either under- 
going exponential growth to much larger masses and higher 
luminosities, or in sub-Eddington quiescent states going into 
or coming out of a period of peak activity. The "break" in the 
QLF corresponds directly to the maximum in the intrinsic dis- 
tribution of peak luminosities, which falls off at both brighter 
and fainter luminosities. The faint-end QLF slope is then de- 
termined by the faint-end slope of the luminosity-dependent 
quasar lifetime (i.e. the differential time the quasar spends in 
any given luminosity interval), for quasars with peak lumi- 
nosity near the observed break. In other words, the number 
of quasars in efficient, near-peak growth stages needed to ac- 
count for the observed bright end of the QLF is already suffi- 
cient to account for the observed faint end, as any given quasar 
is much more likely to be observed at luminosities well be- 
low its peak, and thus the probability of seeing such brighter 
sources at these lower luminosities entirely determines the 
shape of the faint-end QLF. 

This modeling of the quasar lifetime, and the corresponding 
interpretation of the QLF, then provides a physical motivation 
for the location of the break luminosity a nd the faint-end slope 
of the o bserved QLF (Hopkins et al. 2005c d e). Furthermore, 
iHopkin s et al. (2005e) demonstrate that this description accu- 
rately predicts a large number of quasar observations, includ- 
ing the QLF in all wavebands and at all redshifts, the dis- 
tribution of Eddington ratios as a function of luminosity and 
redshift, the BH mass function (both of active and relic BHs), 
the distribution of column densities and obscured fraction as 
a function of luminosity, the X-ray background spectrum, and 



the anti-hierarchical BH growth described above. Given this 
remarkable agreement and the simple physical motivation for 
the faint-end QLF slope unique to this modeling, it should 
account for both the value and evolution of this slope with 
redshift. 

In this paper, we use our picture of quasar lifetimes and the 
resulting QLF interpretation from Hopkins et al. (2005a-f) to 
predict the faint-end slope of the quasar lifetime as a func- 
tion of quasar peak luminosity, which in this model can be 
immediately translated to a prediction of the faint-end QLF 
slope as a function of the QLF break luminosity and, conse- 
quently, as a function of redshift. We are able to accurately 
predict the observed evolution in the slope and corresponding 
LDDE based only on the input physics of our hydrodynamic 
simulations, as a unique consequence of our interpretation of 
the QLF and without the need to invoke any evolution in the 
shape of e.g. the quasar host or galaxy mass distribution. We 
further develop a simple analytical model of quasar feedback 
which gives similar predictions, and demonstrates the scale- 
dependent feedback physics which drive the evolution of the 
faint-end slope. 

Throughout, we define L = Lboh the bolometric luminosity. 
We denote the faint-end slope 7 by = d^/dlog{L) oc L~^, 
where 4> is the differential bolometric QLF and we are con- 
sidering luminosities L L^,, where is the "break" in the 
QLF. Note that d^/dL cx L~^~\ as e.g. in the standard "dou- 
ble power law" form of the QLF. The sign convention (i.e. 
positive 7) also follows the general observational convention 
and standard double power law convention, with d^/dL cx 
1 /{(L/Lo)'"^' + (L/LoV-^^ }. We adopt a I^m = 0.3, Qa = 0.7, 
Hq = 70kms"' Mpc"' cosmology. 

2. THE FAINT-END SLOPE AS A FUNCTION OF PEAK LUMINOSITY 

If quasars spend a differential time dtq/ d\og{L) per loga- 
rithmic interval in luminosity, then the observed QLF is (at 
times when the quasar lifetime is short compared to the Hub- 
ble time) 

<j){L) = J dtQ/d\og{L)ri{Lps^i,)d\ogL (1) 

where n(Lpeak) is the rate at which quasars of a given peak 
luminosity Lpeak are created or activated (per unit time per 
unit comoving volume per logarithmic interval in Lpeak)- In 
Hopkins et al. ( 2005(® , we use this and the realistic model- 
ing of quasar lifetimes from 'Hopki ns et afl J2005alhleh to de- 
termine n(Lpeak), and find that it does not have the same shape 
as the observed QLF, as is predicted by simple models in 
which the quasar turns "on/off" as a step function or follows 
a pure exponential light curve. Instead, «(Lpeak) traces the 
shape of the observed QLF at the bright end (above L*), peaks 
at ipeak ^ L^,, then falls off below this. Therefore, the slope 
of the faint-end QLF, 7, is dominated by the faint-end slope 
of dtq/ d\og(L) for quasars with Lpeak ^ L*; i-e. the observed 
faint end of the QLF is dominated by brighter peak luminosity 
sources, accreting efficiently but in early merger stages or in 
sub-Eddington states transitioning into or out of peak quasar 
activity, so 7 is determined by the integrated probability of 
seeing the population of brighter Lpeak sources, dominated by 
sources with Lpeak L*, at lower L. 

Given this interpretation of the QLF, 7 should be simply 
calculable from the L -C Lpeak slope of dtq/ d\og{L) (since 
the probability of seeing such a source at L is proportional to 
its lifetime at L), for sources with Lpeak ~ L* . If the quasar 
lifetime is a function of L and Lpeak only (i.e. not affected 
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systematically by other host galaxy properties), then 7 can 
be entirely predicted by knowing which directly gives 

the peak in the «(Lpeak) distribution as a function of redshift. 
There may be some curvature introduced because of the non- 
zero contributions to the faint-end of the QLF from sources 
with Lpeak T^i*, but as such corrections depend on the exact 
shape of «(Lpeak) and are small, as n(Lpeak) drops off rapidly 
(as the bright-end QLF does) away from L = L*, we ignore 
them here. 

2.1. The Quasar Lifetime from Simulations 

We first consider the L ^ Lpeak behavior of the quasar 
lifetime determined from hydrodynamical simulations. We 
consider a suite of sev eral hundred simulations, described in 
detail in iRobertson et al. (2005) and Hopkins et al. (2005e). 
These were performed with the new parallel TreeSPH 
code GADGET-2 (Springel 2005), which uses an entropy- 
conserving formulation of smoothed particle hydrodynam- 
ics (SPH; S pringel & Hernquislll2002l) . and includes a sub- 
resolution, multiphase model of the den se interstellar medium 
(ISM) to describe star formation (Springel & HernauisJ 
12003). The multiphase gas is pressurized by feedback from 
supernovae, allowing us to stably evolve even pure gas disks 
(see, e.g. Springel et al. 2005b, Robertson et al. 2004< 
Springel & Hernquist 2005). BHs are represented by "sink" 
particles that accrete gas, with an accretion rate M estimated 
using a Bondi-Hoyle-Lyttleton parameteriza tion, with an up - 
per limit equal to the Eddington rate ( Spri ngel et alJl20'05bh . 
The bolometric luminosity of the BH is then L = erMc^, where 
e = 0.1 is the radiative efficiency. We further allow a small 
fraction (typically w 5%) of L to couple dynamically to the 
gas as th ermal energy. This frac tion is a free parameter, deter- 
mined in lDi Matteo et al J (120051) by fitting the Mbh-ct relation. 
We generate two stable, isolated disk galaxies, each with an 
extended dark matter halo having a Hernquist ( 1990) profile, 
an exponential disk, and a bulge. The galaxies are then set to 
collide from a zero energy orbit. 

In our suite of simulations, we vary the masses and virial 
velocities of the initial galaxies, halo concentrations, ISM 
equation of state, parameters describing feedback from super- 
novae and black hole growth, presence or absence of bulges 
in the host galaxies, initial BH seed masses, numerical resolu- 
tion (typically ~ 2 x 10^ particles per galaxy, but we consider 
up to 10-100 times as many), disk inclinations and pericenter 
separation of the initial orbits, and initial disk gas fractions. 
We further scale all galaxy properties appropriately to resem- 
ble galaxies at redshifts z = - 6, for a large sub set of our 
simulations as described in Robert son et alJ (120051) . Our sim- 
ulations produce quasars with L, Lpeak from ~ 10^- IO'^'Lq, 
spanning the entire range o f observ ed quasar luminosities at 
all redshifts. iHopkins et alJ (|2005e) consider this set of sim- 
ulations in detail, and use it to determine the quasar life- 
time dtg/dlogiL). Critically, they find that dtQ/ d\og(L) 
shows no systematic dependence and little scatter with any 
of the varied parameters, whe n quantified as a fun ction of L 
and Lpeak- The predictions in Hop kins et alJ (12005 d) did not 
depend sensitively on the (much more uncertain) faint-end 
slope, and therefore the authors did not consider a detailed 
fit to the faint-end behavior of dtq / d log(L), but rather param- 
eterized the lifetime as a simple exponential, dtg / d log(L) = 
fo exp(-L/Lo), where Lq ~ 0.2Lpeak and fo depends weakly on 
Lpeak, which providcs an excellent fit to the simulation results 
for L > 10"2- 10"'Lpeak. iHopkins et all (l2005bh consider the 



L <C Lpeak behavior of dtq/ d\og(L) in more detail, and find 
that the lifetime more closely resembles a power law at these 
L, with a power-law slope a a function of Lpeak- The combi- 
nation of these results suggests that the quasar lifetime is best 
parameterized as a Schechter function with slope a, normal- 
ization fo(Lpeak), and turnover Lo(Lpeak), for purposes where 
the faint-end slope is important; i.e. (to clarify our conven- 
tions) 

dlog(L) "VLo/ L,,/ 

For a QLF with break L*(z), which implies an n(Lpeak) distri- 
bution peaked at Lpeak ^ L^, the observed faint-end QLF slope 
is then7«Q;(Lpeak = L*). 

Here, we consider the faint-end (L ^ Lpeak) slope from 
these simulations in greater detail. For each simulation, we fit 
a Schechter function to the quasar lifetime as a function of lu- 
minosity. We find similar results to Hopkins et al. (2005e) for 
Lo and fo (allowing for the well-known degeneracy between 
fo and a in these fits), but n ow quan t ify a as a function of 
Lpeak- As discussed in IHopkins et alJ (|200531), these faint- 
end slopes are subject to several uncertainties in our mod- 
eling (the reason for ignoring the correction at very low lu- 
minosities implied by a in Hopkins et al. 2005e). One is 
the finite time duration of our simulations, which may flat- 
ten dtq/ d\og{L) at low L, as the BH cannot completely relax 
in the given time. A complementary means to determine 7 is 
then to consider the rate at which L falls off after L = Lpeak 
at time t = tp = f (Lpeak)- We take L(f) /Lpeak to be a function 
of t-tp, and find that it falls off in approximate power-law 
fashion for all our simulations. Lit) oc (t-tp)~^. This implies 
dtg/ d\og{L) cx L~^l^ , i.e. 7 = 0;=! / [3. Although these fits are 
not affected by the finite time duration of the simulations, they 
do not include the time spent at different L before L = Lpeak, 
and thus are entirely accurate only at f » tp or for a symmetric 
rise/fall in L, and may therefore overestimate the steepness of 
7, so we consider both fitting methods below. 

The uncertainties in this mo deling at low luminosi - 
ties are discussed in detail by IHopkins et alJ (I2005cleft . 
IHopkins et ai] ( l2005bh examine fits to the power-law behavior 
of dtg/ d\og{L) at L ^ Lpeak considering both the entire sim- 
ulation and only times after the merger (similar to our fits for 
the decay of L(f)). Further, H opkins et al. ( 2005 c) consider 
the application of an ADAF-type correction for radiatively in- 
efficient accretion flows at low accretion rates (following e.g., 
Naravan & Yi 1995), to account for more detailed variation 
in the radiative efficiency and spectrum as a function of ac- 
cretion rate. We apply this correction in our fitting as well. 
Although these previous works consider a much smaller sub- 
set of the simulations we have fit to, they find similar results, 
and we compare their fits and ours below to demonstrate both 
the general agreement and range of uncertainty introduced by 
these different determinations of 7. 

Figure ^ shows the results of this fitting for several rep- 
resentative simulations. We show the lifetime as a func- 
tion of L (for clarity, the integrated lifetime above a given 
L, tQ{L' > L), is shown), for simulations with Lpeak /Lq = 
3.1 X 10'^ (black circles), 7.3 x 10'^ (blue diamonds), 1.7 x 
10'^ (green squares), 1.4 x 10'"^ (red x's), and the best-fit 
Schechter function (solid lines) to each, with slope a, on the 
left. We compare these to the slopes 1/(3 determined by fit- 
ting L(f)/Lpeak as a function of t-tp (right, same notation). 
The resulting faint-end slope 7 for each Lpeak and each fitting 
method are also shown; there is a clear decrease of 7 with 
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Fig. 1 . — Left: Quasar lifetime (integrated time above a given L), for representative mergers with peak luminosities Lpeak/i© = 3.1 X 10'^ (black circles), 7.3 X 
10'- (blue diamonds), 1.7 X lO" ( green squares), 1.4 X lO''' (red x' s), with the best-lit Shechter functions (lines) of slope a. Right: Bolometric luminosity 
over peak luminosity (L/tpcak) as a function of time (tf = f-f(Lpcak)). with the best-fit power law of slope /3. The predicted faint-end LF slopes 7 = Q = l//3are 
also shown. 



TABLE 1 
Faint-End Slope 7(Lpeak) 



Determination" 


7l2=7(Vak = lO'2L0) 


rflogCipciik) 


H05b, all times 


0.95 ±0.25 


-0.32 ±0.11 


H05b,f >f(Lp,^) 


0.76 ±0.20 


-0.30 ±0.07 


H05c, all times 


0.94 ± 0.20 


-0.33 ±0.08 


H05c,r >f(Lpeak) 


0.50±0.15 


-0.50 ±0.15 


Schechter fitting 


0.44 ±0.02 


-0.21 ±0.02 


L(t) oc decay 


0.69 ±0.05 


-0.30 ±0.04 


Cumulative best-fit 


0.55±0.12 


-0.25 ±0.04 


Blastwave model 


0.35 ±0.15 


-0.22 ±0.06 



"H05b, H05c (Hopkins et al. 2005b,c) use a small subset of the simulations 
here, with H05c applying an ADAF-type correction at low accretion rates. 



^peak • 

We fit for 7(Lpeak) in each of our suite of several hundred 
simulations with both methods described above. In all cases 
we find 7 decreases with Lpeak in a manner similar to that in 
Figure Given the scatter, the resulting slopes 7(Lpeak) are 
well fitted by a simple log-linear relation, 

7(^peak) = ln+-p. 77 T (log(Lpeak/i0)- 12), (3) 

<ilog(Lpeak) 

where 712 ~ 0.5 and -n-rr — r ~ -0.3. This is illustrated 

"■^ dlog(Lp^.ak) 

in Figure |2] where the value of 7 determined from fitting 
Schechter functions to the quasar lifetime for each of our sim- 



ulations is plotted as a function of the peak luminosity of those 
simulations (circles, with errors), with the best-fit log-linear 
relation (line) (see also Figure 7 of Hopkins et al. 2005b). The 
results for both fitting met hods are summarized in TablelD We 
also show the results of Hookins et al. ' ('2005b"c"), described 
previously. As shown in detail in Hopkins et al. (2005e), we 
find no systematic dependence on redshift, host galaxy prop- 
erties, or any other varied parameters in the simulations. 

In our simulations, gas is expelled when feedback from ac- 
cretion heats the surrounding gas to the point where it can 
no longer efficiently cool and rapidly unbinds it. Around 
more massive BHs (i.e. higher-Lpeak quasars), this "blowout" 
event is progressively more violent; the higher luminosity of 
the quasar heats the gas more rapidly and to higher tempera- 
tures, resulting in a more violent and complete expulsion of 
gas (see also Cox et al. 2005, who show the gas mass fraction 
expelled increases with BH mass). This is clear in the rep- 
resentative simulations for which we fit L(t-tp), which falls 
off much more rapidly in higher-Lpeak simulations as the den- 
sity decreases and temperature (sound speed) increase more 
dramatically around the quasar 

2.2. A Simple Model of Quasar Feedback 

To understand this trend in greater detail, we consider an 
analytical model for the quasar "blowout" event, i.e. the rapid 
heating and expulsion of gas when the quasar reaches a criti- 
cal mass/luminosity, in agreement with the Mbh-ct relation. 
We can use such a model to predict and understand at least the 
late-time (post critical mass) evolution of the quasar accretion 
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logioC-peak/L®) 

Fig. 2. — Faint-end QLF/quasar lifetime slope 7 = determined from fits to the quasar lifetime as a function of peak luminosity Lpeak in our simulations. The 
best-fit log-linear relation (see TablelT] is shown (solid line). 



rate, although the early activity (which as demonstrated above 
is a significant contributor to the lifetime at low-L) is deter- 
mined by the chaotic merger process and is much more diffi- 
cult to accurately predict in simple analytical form. We note 
that this problem has been considered in greater detail in pre- 
vious works (e.g. Barkana & Loeb 2001; Furlanetto & Loeb 
|2p01; Scannapieco & Oh 2004), but these authors have con- 
sidered the large-scale structure and impact of these outflows 
on the outer regions of halos and clusters, the surrounding 
IGM, and reionization; we are instead interested in the very 
small scale, rapid expulsion of gas from the innermost BH- 
dominated accretion regions of the galaxy. Interestingly, these 
works suggest that feedback on these scales is critical in regu- 
lating the slope of the bright end of the QLF, which highlights 
the importance of quasar feedback on all scales in the QLF. 

In our numerical modeling ( Di Matteo et al. 2005 ) and most 
analytical models of the Mbh-o' relation (e .g.. Silk & Rees 
pQS ' Fabian 1999: iCiotti & Ostrikedl2()()U .Wvithe & Loeb 
l2002h . feedback energy from accretion onto the central BH 
heats the surrounding gas (or momentum from coupling of 
the radiation field to dust drives a wind). In earlier stages 
of BH growth, this is a small luminosity, and the gas is able 
to efficiently cool and re-radiate this energy. However, in the 
peak merger stages, the BH growth is exponential (Eddington- 
limited), and as it rapidly grows in luminosity, a critical 
threshold is reached where the nearby gas can no longer cool 
efficiently in a local dynamical time. As a result, the gas is 
rapidly expelled, and what remains is heated to virial temper- 
atures, preventing it from being easily accreted, and further 
accretion shuts down. The threshold for this behavior is deter- 



mined by the local gas properties and gravitational potential 
of the host, but essentially all models of feedback driven self- 
regulation (owing to exceeding a critical energy input in of 
order a dynamical time) predict a relation between final black 
hole mass and virial velocity or spheroid velocity dispersion 
of roughly Mbh oc V^^^. or Mbh oc a*, respectively. 

Because the gas in the central galaxy regions is able to effi- 
ciently cool prior to this self -regulation, and is then suddenly 
heated and driven out in a short time (in our simulations, this 
occurs over a timescale ^ 10^ yr in massive mergers) by an 
exponentially increasing luminosity, it is reasonable to model 
this dri ven outflow phase as a Sedov-Taylor-type blastwave 
JSedovt 1 946- . 1 959; Tavlor..l 950i) . with energy injection from 
a point explosion with energy E = riLLEddiM^u^ t^yn {tjl de- 
scribes the efficiency of feedback coupling, « 5% in our sim- 
ulations and similar in most analytical models of the Mbh-c 
rela tion). A more detailed e xamination of this problem in 
e.g. lFurlanetto & LoebH2001h shows that this is a remarkably 
accurate approximation to a full solution including radiative 
cooling, the pressure of the external medium, magnetic fields, 
and further effects. Modeling the accretion rate with a Bondi- 
Hoyle (Bondi&Hovle 1944; Bondi 1952) parameterization 
(as by definition we are considering times when the accretion 
rate falls below the Eddington limit), considered at the radius 
of influence of the BH, Rbh = GMbh<^~^ (c being the spheroid 
velocity dispersion) gives L oc M|jjp(7?bh)c^''(/?bh), where 
is the isothermal sound speed (in our simulations, the effec- 
tive sound speed of the multi-phase ISM is used, see Springel 
& Hernquist 2003). As the accretion rate is rapidly falling 
from its peak during this stage, we can reasonably approxi- 



6 



Hopkins et al. 



mate Mbh = -^bh =constant, and thus /?bh =constant, and we 
only need to calculate the time evolution of the gas density 
p and sound speed c, at fixed r = 7?bh to determine the time 
evolution of luminosity L in this stage. We are only interested 
in the scaling of L{t) oc (where t here is the time since the 
peak quasar luminosity or beginning of this "blowout" stage), 
so we could theoretically pick any fixed radius to evaluate 
these quantities, and thus the choice of the radius of influence 
of the BH as opposed to the transonic radius, for example, is 
not significant. 

We first consider the simplest, scale-invariant solution to 
this problem, in which we neglect the gravitational fi eld of 
the BH and halo, following lOstriker&M^K^ lTmS^ m our 
derivation. We also ignore the consequences of radiation as 
noted above, but this should not be a large effect in the very 
inner regions of the blastwave expansion in which we are in- 
terested (see also Murray et al. 2005, who consider in detail 
the coupling of BH radiation to dust and gas opacity and show 
that it produces a very similar Mbh - c relation and quali- 
tative "blowout" behavior). As we need to describe p{r,t) 
at small radii (relative to the galaxy size and final radius of 
the blastwave), we must consider the density gradient of the 
medium, for which we adopt the one-power approximation 
(OPA) p oc r"*^*" for the ambient (pre-shocked) medium. There 
may be some non-spherical effects due to asymmetry in the 
density profile or e.g. bipolar outflows driven by the quasar, 
but the result is qualitatively quite similar to the spherical so- 
lution along any given sightline, and we furthermore expect 
that the chaotic interactions and short dynamical times in the 
central regions of the merger will rapidly isotropise the region 
of interest. For a self-similar OPA Sedov-Taylor solution, the 
shock radius expands as Rs oc f with the post-shock (in- 
ternal) density profile p = pi (r/Rs)''', where pi oc p{Rs) is the 
density just inside the shock front. At fixed r = Rbh, then, we 
have 

p(i?BH,f) (4) 

For a self-similar, energy-conserving blastwave, rj = 2/(5 - 
kp), and mass conservation (coupled with the strong shock 
jump conditions;) requires Ip = (6-(7eos+ l)^p)/(7E0S- !)■ 
Here, we denote the ratio of specific heats as 7eos to distin- 
guish it from the faint-end QLF slope 7. We then obtain 



(kp + lp)r]-- 



■l\5-kJ 



7eos- 



(5) 



The pressure P = pel also follows a power-law scaling (with 
internal power-law slope Ip and external slope kp), but with a 
more complicated result for Ip, 

37ios + 2O7EOS + 1 - (7E0S + 1 )(37eos + 1 



lp = 



2(7ios-l) 



(6) 



Thus the final scaling of L{t) oc pcj' oc p^l^P in this phase 
is given by L oc f with 



7eos- 



(7) 



This formalism has an exact, self-consistent solution corre- 
sponding to e.g. a blastwave in an isothermal sphere or wind, 
which should be a reasonable approximation to the mergers 
our simulations model in detail, with kp = 2 and 7eos = 5/3 
(yielding Ip = \, Ip = 3). This solution gives (3 = 2, imply- 
ing a faint-end QLF slope 7 = 1 /2. But the above equati ons 
should be a good approximation JOstriker & McKeell988h for 



a range of 7eos ^ 1-5/3 (corresponding to variability in the 
equation of state owing to the density and temperature depen- 
dence of star formation and radiative cooling) and kp ^ 1-3 
(corresponding to e.g. the inner regions of a Hernquist 1990 
spheroid profile with a central point mass). In fact, it is still 
exact for the equilibrium structure of a filled blastwave with 
the most natural v oc r velocity structure, which yields 



kp — kp crit ~ 



7 - 7eos 
7E0S + 1 ' 



f3 = 



1 



37EOS-1 



20--(7Eos+l) 



(8) 



(9) 



For the reasonable range of equations of state 1 < 7eos < 5/3, 
this results m2 < (3 < 11/2, i.e. 0.2 ;< 7 < 0.5, in good agree- 
ment with the range of 7 observed and determined in our sim- 
ulations. 

The scale-invariant Sedov-Taylor solution then explains 
well the typical values of 7 implied by our simulations and 
the rapid falloff of L{t) during the "blowout" phase of quasar 
evolution. However, this does not immediately describe the 
dependence of 7 on Lpeak- This dependence is not contained, 
as might be presumed from the above derivation, in 7eos, as 
essentially all our simulations have 7eos in the entire range 
^ 1-5/3 in their central regions and we specifically find 
no significant systematic dependence of 7 or the quasar life- 
time on the ISM gas equation of state (see also Hopkins et al. 
2005e, Figures 4 & 5). Rather, the dependence of 7 on Lpeak is 
driven, naturally, by the fact that this problem is not precisely 
scale invariant or self-similar 

To illustrate the dependence on Lpeak, consider a small 
change in the inital logarithmic slope of the density in the 
central regions, kp = k^p + 5kp. For simplicity we expand about 

the k^p = 2, 7eos = 5/3 exact solution, but our result is similar 
regardless of the choice of these parameters (within the ranges 
defined above). This gives 



l3 = 2-n5kp/3, 



(10) 



(11) 



i.e. d^/dkp = 13/12. Thus, if the initial logarithmic density 
gradient in the inner regions is steeper, i.e. has a larger kp 
(po oc r~'^>'), then (3 is actually shallower. In detail, two com- 
peting effects occur. First, the density falls off more rapidly, 
allowing the blastwave front to propagate more rapidly as it 
builds up more of its mass earlier and encounters less mass 
further from the central regions, but this effect is relatively 
weak, i.e. ri + 26kp/9. The dominant effect is the alter- 
ation of the post-shock density profile. Because the density 
gradient is steeper, the propagation of the blastwave builds up 
a less pronounced "snowplow", i.e. the mass buildup at the 
wavefront is less pronounced, implying a flatter post-shock 
density gradient by mass conservation (recall that p increases 
with radius within the shocked region). Essentially, with the 
gas mass concentrated in the center, the shocked density pro- 
file evolves less dramatically, because less mass is added and 
the early blast is able to effectively redistribute more of the 
ultimately acquired mass. Given that less massive halos typ- 
ically have higher concentrations, we can estimate the mag- 
nitude of this effect on the evolution of 7 with Lpeak- Esti- 
mating dkp/dc ^0.1 (where c is the concentration index and 
we estimate a change 10 in c changes the inner logarithmic 
slope at fixed r by ~ 1), we can determine Jc/ii log(Lpeak) 
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using Lpeak Mbh oc Mvir (e.g., Marconi & Hunt 2003) and 
c « 9(Mvi,/lO'3M0)-o 13 (iBuUocket al.l2001l) . and we find 

-Jl 2f-^y°" (12) 

rflog(Lpeak) VlO'ZLeJ ' ^'^^ 

in reasonable agreement with our measured dependence 
of 7(Lpeak) in our simulations and only weakly dependent 
on Lpeak- This actually predicts that the magnitude of 
d7/JlogLpeak should be larger at low-Lpeak ('^ -0.28 at 
ipeak --10" Lq) and smaller at high Lpeak, (^ -0.16 at Lpeak ^ 
IO'-'Lq) which may occur (see Figure |2j, but this is most 
likely a coincidence, as our modeling of the blowout in scale 
invariant fashion and the application of a concentration pa- 
rameter in such chaotic period in evolution of the merger are 
rough approximations at best. 

In this model of the quasar blowout phase, however, the ap- 
proximate values for the rate at which accretion declines are 
simply determined by a standard Sedov-Taylor solution in a 
medium with a density gradient. The sudden injection of en- 
ergy as the BH crosses a critical mass/luminosity threshold 
and the surrounding gas is no longer able to efficiently cool 
drives a strong outflow and heats the remaining gas, rapidly 
shutting down accretion. Incorporating the weak, but not 
negligible, effects of the change in concentration with mass 
breaks the scale invariance of this solution, as lower-mass sys- 
tems have steeper inner density profiles, which flatten the evo- 
lution of the accretion rate in time and produce steeper 7. Of 
course, many other effects will break the self-similarity of this 
problem as well - a realistic gravitational potential will imply 
a characteristic scale length, and the physics of radiative cool- 
ing will likewise define fundamental physical scales. Even 
the scale-invariant solution incorporating radiative energy loss 
depends on the logarithmic slope of the cooling function vs. 
density and temperature, but the values of these slopes depend 
themselves on the characteristic temperature of the blastwaves 
and change quite significantly over the mass scale of our sim- 
ulations (heating to virial temperatures implies tem- 
peratures T ^ 10^- 10^ K over the mass range shown in e.g. 
Figure|2ll. Although we do not model the chaotic interactions 
and evolving BH mass of the early merger stages, the more vi- 
olent torquing associated with more massive mergers can ex- 
plain the more rapid, peaked BH evolution over a larger range 
in BH mass, even in early merger stages, generating a flatter 
quasar lifetime which spans a wider range in luminosity. 

Given these various scalings, it is possible that our observed 
trend of 7(Lpeak) could change or even reverse at very low 
Lpeak, as in models of stellar winds in dwarf ellipticals for 
which lower-mass objects (lower Mbh) are more easily un- 
bound (e.g. Mac Low & Ferrara 1999) (although this is more 
concerned with the large-scale binding of gas, as opposed to 
evolution in the inner accretion regions of interest in our mod- 
eling), but the masses/luminosities where this is likely to be- 
come important (Mgai < 10** Mq, i.e. Lpeak ^ 10^ L0) are well 
below the break luminosity at any redshift, and thus will not 
affect our results. Likewise, this could occur at very large ra- 
dius r a in any Lpeak system, but again we are not attempt- 
ing to model the large-scale blastwave but only the evolution 
relevant to the small accretion regions. 

2.3. Parameterizations of the Quasar Light Curve 

The fits and analytical modeling above imply a useful, sim- 
ple prescription for the quasar light curve as applied in semi- 
analytical models and other theoretical modeling which can- 
not resolve the detailed time history of individual objects as 



we can in our simulations. Generally, the quasar light curve 
is characterized by two "modes": a "growing mode," charac- 
terized by high-Eddington ratio, rapid black hole growth, and 
a "decaying mode," characterized by the nearly self-similar 
power-law falloff of the quasar luminosity as nearby gas is 
heated or expelled. The "growing mode" can be most simply 
parameterized by exponential growth at a constant Eddington 
ratio m, with an exponential light curve L = L{t = 0) exp {t / tg], 
where tg = ts/m is the e-folding time and ts = 4.2 x 10^ yr is 
the Salpeter time. Such black hole growth is expected in es- 
sentially all models of quasar activity, where plentiful fueling 
easily enables high accretion rates. 

Once the quasar reaches a critical luminosity or mass (in 
e.g. accord with the M-a relation), it begins to heat and ex- 
pel the surrounding gas, and the accretion rate rapidly falls off 
- i.e. the light curve can be roughly parameterized as enter- 
ing the "decay mode" described above. Our fits to the quasar 
lightcurves after the "blowout" phase and analytical model- 
ing of this phase of evolution as a driven blastwave suggest 
a power-law decline in the quasar lightcurve, which can be 
simple modeled as 



L(f)/i'peak ~m 



1 



(13) 



Lpeak here is the peak bolometric luminosity, just as the quasar 
enters the "blowout" phase, and L(f ) is the subsequent bolo- 
metric luminosity at time t after the peak (L = Lpeak at f = 0). It 
is also not a bad approximation to use this to model both the 
light curve and accretion rate with L(f)/Lpeak ~ because 
not much total black hole mass is accumulated in this mode. 
The equation shown assumes m « 1 at peak luminosity, as is 
true in most of our simulations, but this can easily be renor- 
malized to any assumed constant m in the constant Eddington 
ratio "growing mode." The functional form of this equation 
is chosen such that it joins continuously with the constant Ed- 
dington ratio exponential light curve at f = (L = Lpeak, i-C- 
at the beginning of the "blowout" stage), and behaves as our 
fitted power-laws at times large compared to the duration of 
the "blowout." 

The critical behavior determined from our simulation is 
this power-law decay of the quasar light curve at late times, 
L(f) oc (t/tg)-^. We have measured (3 directly in our simu- 
lations above, and estimated its value from simple analytical 
models of the quasar-driven "blowout" phase. To lowest or- 
der, a canonical value of /3 = 2 is suggested by both our fits to 
the simulation light curves and the self-similar Sedov-Taylor 
solution for a quasar-driven blastwave and Bondi accretion. 
But we have also explicitly determined /3 as a function of peak 
luminosity, which we inverted to determine 7(Lpeak) above. 
Our fits to the simulations yield 



d(3 



C/ log Lpeak 



{Iogio(ipeak/1O'2L0)}. (14) 



where 1.7±0.3 and t//3/t/logLpeak « 0.7±0.1. There 

is some ambiguity in (3 depending on whether we are fitting 
to a given blowout or a whole quasar lifetime which includes 
the "growing" phases (giving shallower effective jS), as can 
be seen directly by the fact that different values of the faint- 
end quasar lifetime slope 7 are suggested by our fitting of 
the quasar lightcurves and our fitting of the quasar lifetime to 
a Schechter function. However, for most theoretical models 
which require a time-dependent quasar lightcurve, the param- 
eterization desired is that of an individual blowout, and so 
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these fits are used to determine Equation above. Further- 
more, the "effective" (3 from fitting to the quasar lifetime is 
simply \/a, which can be directly determined from the values 
of 7(Lpeak) given in Tabled There is a significant amount of 
scatter around this mean relation in (3 (in contrast to the quite 
small scatter in a) at a given Lpeak, A/3 ~ 1, which can easily 
be incorporated in theoretical models of quasar populations if 
desired. Clearly, this formula for /3(Lpeak) cannot be extrapo- 
lated to arbitrary luminosities, as by definition /3 > always. 
In fact, most of the simulations which are well-fit by a very 
low (3 are so because they never reach a high Eddington ratio 
and cause a significant "blowout" event, but rather they quies- 
cently accrete for ~Gyr at moderate (~ 0. 1) Eddington ratios. 
In either case, this is only imporant Lpeak ^ a few x IO^Lq, 
which implies black hole masses < IO^Mq, well below the 
regimes where the processes we have modeled should be im- 
portant (and below the limits to which our simulation light 
curves can be reliably estimated). 

There is, unfortunately, significant ambiguity in the appro- 
priate value of fg to use in this simple parameterization of the 
quasar light curve. If we use Equation^Jto estimate the cu- 
mulative quasar lifetime at low luminosities (L <C Lpeak, low 
enough that the contribution from the "growing mode" is neg- 
ligible) and compare this with the Schechter function fits or 
direct calculation of the quasar lifetime from our simulations 
at these luminosities, there is a well-defined approximately 
constant f g « f ^ = 4.2 x 10^ yr (or, if we allow fg to vary with 
P and demand that this match to the low-luminosity limit of 
the fitted Schechter function with the same Lpeak, tg « (3ts ^ 
10** yr). However, there are several minor "blowout" events in 
our simulations, and usually at least two major ones (one at 
first passage, the other after the final merger). The fg which 
we have constrained by fitting the integrated quasar lifetime 
is approximately the sum of these contributions. If we instead 
fit to fg in an individual blowout, we obtain a significantly 
smaller value fg w 10^ yr (as is also suggested by direct ex- 
amination of Figure^. The "appropriate" value of fg is then 
somewhat determined by the ability of a given model to re- 
solve separate starbursts triggered by e.g. passage of galaxies 
as opposed to the final merger, and resolution of triggering 
by minor mergers. However, our simulations constrain the 
range of reasonable fg to ^ 10^- 10** yr, which is still a sig- 
nificantly more narrow range than the observations constrain 
(e.g., ^Martini 2004). More important, to first order, fg only 
controls the normalization of various predicted quantities and 
it should be easily calibrated in a given theoretical model by 
comparison to e.g. quasar number counts and the normaliza- 
tion of the luminosity function. Even with these uncertainties 
in the proper "effective" fg due to the limited resolution of 
individual events in semi-analytical or cosmological models, 
the implied appropriate effective fg for a given model is an 
interesting constraint on e.g. radiative efficiencies of quasars 
(as the fg we fit scales with the e-folding time for black hole 
growth, and hence the radiative efficiency) and the effective- 
ness of quasar feedback. For example, a very short fg, es- 
pecially in a model which is effectively summing over sev- 
eral different "blowout" events, implies very rapid expulsion 
or heating of gas and therefore highly efficient coupling of 
quasar energy or momentum to the surrounding gas in the 
"blowout" phase. 

The quasar lightcurve is, of course, much more complex 
than we have modeled here. Not only are there several 
such "growing," "blowout," and "decaying" modes in a given 



merger event, but each follows a light curve which is not triv- 
ial in detail, and our mixed exponential and power-law light 
curves are not always a good fit. Whenever possible, theoret- 
ical models should adopt the most accurate approximations 
to the quasar light curve available, attempting to resolve the 
detailed time history where it is important (for example, in 
estimating the properties and luminosity function of the faint 
quasar population). However, in many cases this is not pos- 
sible or feasible, and a simple analytical approximation for 
the quasar light curve is necessary or convenient. The above 
simple parameterization describes both the growing and de- 
caying phases of the quasar light curve with reasonable accu- 
racy, and captures the critical properties of the quasar lifetime 
and its increase with decreasing luminosity and the tendency 
of bright peak luminosity quasars to spend extended periods 
in low-luminosity phases, and should therefore be useful in 
such modeling. This does not, however, include the effects 
of obscuration, which are not important (at least along most 
sightlines) in e.g. the late blowout stages, but can dramatically 
change the observed quasar light curve in various bands near 
peak luminosity, especially before the "blowout" phase (e.g., 
Hopkins et al. 2005b e). These effects can be estimated from 
the parameterization of column density as a function of im- 
mediate and peak luminosity in Hopkins et al. (2005e), but as 
noted therein, these fits are only approximate and do not accu- 
rately describe the "blowout" stage in which luminosity and 
column density both vary rapidly. 

3. DISCUSSION & CONCLUSIONS 

Realistic, luminosity-dependent quasar lifetimes, in which 
quasars spend more time at lower L, imply a novel interpre- 
tation of the QLF (Hopkins et al. 2005c), in which the faint- 
end is composed of sources with larger Lpeak in early growth 
or sub-Eddington states. The faint-end slope 7 is then deter- 
mined by the slope of the quasar lifetime vs. L, for quasars 
with Lpeak ^ L^, the break in the observed QLF correspond- 
ing directly to the peak in the distribution of quasar peak lumi- 
nosities (BH masses) being formed at any redshift. We have 
determined 7 from both a large range of hydrodynamical sim- 
ulations (varying BH and host galaxy properties) and from a 
simple analytical model of quasar feedback and the Mbh-c 
relation, and find that it is a monotonic decreasing function of 
Lpeak over the range of observed L*. 

As the break luminosity L*(z) is well-known from observa- 
tions of the QLF, we use our modeling to predict 7(2). Fig- 
ure |3] shows our result for 7 in the redshift interval z = 0-4. 
Above z^2-3, the location of L» becomes uncertain, and we 
consider both pure density evolution (L» = constant or increas- 
ing, solid) and pure peak luminosity evolution (declining L*, 
dashed) at higher redshifts. The yellow shaded range shows 
the range predicted by different fitting methods from Tabled 
lines the cumu l ative bes t-fit. The form of L,(z) is taken from 
[Hopki ns et alJ J2005el f>). but is based directly on fits to the 
observed Q LFs of [U eda et al. ( 12003h . iRichards et a l. (2005), 
and Hasing er et alJ We compare this to observations 

in the hard X-ray (Ueda et al. 2003, blue squares), soft X-ray 
(Hasinger et al. 2005, red circles), optical (cyan; Wolf et al. 
2003, diamonds; Hunt et al. 2004, star; Richards et al. 2005, 
x; Pei 1995, +), and radio (green; Sadler et al. 2002, *; Cira- 
suolo et al. 2005, triangle). We convert these to bolometric 
luminosities and re-fit or rescale 7 with the bolometric cor- 
rections of Marconi et al. (2004), which are also discussed in 
det ail in Hopkins e t al.. ( 2005d e) and are similar to those in 
e.g. iRichards et al.i(l2005l) . Although the uncertainties in the 
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Fig. 3. — Predicted faint-end QLF slope 7 as a function of redshift z, from the fits in Table fll Shaded yellow shows the ^ la range depending on fitting 
method, lines the mean prediction. Above z ~ 2, sohd line assumes pure density evolution for the QLF, dashed line assumes pure peak luminosity evolution. 
Observations are shown (converted to bolometric QLF slopes) from different wavelengths as noted. 



observed faint-end slope 7 are large, we reproduce its value at 
all redshifts. 

Figurel^demonstrates how our prediction for the evolution 
of 7 with redshift translates to a "luminosity-dependent den- 
sity evolution" (LDDE). In each panel, the integrated number 
density <& (per comoving volume) of quasars in each of sev- 
eral luminosity intervals is plotted as a function of redshift 
for z < 5. The upper left shows $ for intervals in bolometric 
luminosity L: logCL/L©) = 9- 10 (black sohd line), 10-11 
(blue dotted), 11-12 (light blue short dashed), 12-13 (green 
dot-dash), 12-14 (orange triple-dot-dash), and 14-15 (red 
long dash). The effects traditionally fitted to LDDE forms 
are clear; the density of higher-L systems rises more rapidly 
to a peak at higher redshift, then falls off more slowly. The 
difference in evolution between two L intervals becomes less 
dramatic with increasing L, as observed. Using the bolomet- 
ric conversions of Marconi et al. (2004) described above, we 
also consider the observed density evolution in other bands. 
The upper right shows density evolution in the B-band in in- 
tervals of B-band magnitude, -20 > Mb > -22.5 (black solid), 
-22.5 >Mb> -25 (blue dotted), -25 > Mb > -27.5 (green 
dashed), -27.5 > Mb > -30 (red dot-dash). This modeling 
immediately demonstrates why observations in the optical 
have not found the dramatic LDDE seen in X-ray samples, 
as the effects do not become pronounced until very low lu- 
minosities not usually probed in B-band samples (as low as 
Mb > -18), and the bolometric corrections actually slightly 
"blur out" the effect as well. The lower left shows our pre- 
diction for the hard X-ray, in three intervals of L(2-10 keV) 
= Lfjx- In order to directly compare with the observations of 



lUeda et all i2003l l). we adopt cgs units and consider the in- 
tervals log(LHx/ergs"') = 41.5-43 (black solid), 43-44.5 
(blue dotted), 44.5-46 (green dashed). For each, we show 
the observed best-fit LDDE model with approximate la er- 
rors, in each observed redshift interval (filled circles of ap- 
propriate color). Likewise, the lower right shows our predic- 
tion for the soft X-ray i n inter vals of L(0.5-2 keV) = Lsx- 
We compare directly to iHasin^er et alj (2005) in the inter- 
vals log(Lsx/ergs"') = 42-43 (black solid), 43-44 (blue 
dotted), 44 - 45 (green dashed), 45 - 46 (yellow dot-dash), 
> 46 (red triple-dot-dash). Again we show the observed best- 
fit LDDE model in each observed z interval (filled circles of 
matching color), but we have multiplied the QLF normaliza- 
tion of Hasinger et al. (2005) by a f actor of 10 to account for 
the mean obscured fraction ( Hopkin s et alJ2005^ . For all the 
above, we adopt pure density evolution for z > 2. We do not 
show the errors in our predictions, as these are actually domi- 
nated by the uncertainty in the fitted n(Lpeak) distribution, not 
the uncertainty in 7(-Lpeak)- 

The results in Figure|4]show that our modeling of the quasar 
lifetime predicts very accurately the "luminosity-dependent 
density evolution" seen in X-ray samples, at all luminosi- 
ties and redshifts. We may slightly underpredict the low- 
luminosity, high-redshift number of sources, but this is unsur- 
prising, as it is both where the observations are most uncertain 
and where, for a fixed low L, the break luminosity becomes 
very much larger than L, meaning that our modeling is based 
on the quasar lifetime at L well below Lpeak where it is most 
uncertain. Still, the agreement is remarkable, and it is de- 
rived almost entirely from the quasar lifetime as a function of 



10 



Hopkins et al. 




012345 012345 
z z 

Fig. 4. — Predicted comoving number density in different luminosity intervals as a function of redshift. Upper left: Bolometric; log(L/L0) = 9—10 (black 
solid line), 10- 11 (blue dotted), 11-12 (light blue short dashed), 12- 13 (green dot-dash), 12- 14 (orange tiiple-dot-dash), 14— 15 (red long dash). Upper 
light: B-band; -20 >Mb> -22.5 (black solid), -22.5 > Mb > -25 (blue dotted), -25 > Mb > -27.5 (green dashed), -27.5 > Mb > -30 (red dot-dash). Lower 
left: Hard X-ray; log(LHA'/ergs"' ) = 41.5 — 43 (black solid), 43-44.5 (blue dotted), 44.5 — 46 (green dashed), compared to observations of Ueda et al. ( 2003) 
(filled circles of corresponding color). Lower right: Soft X-ray; log(Lsx / erg s"' ) = 42 - 43 (black solid), 43 — 44 (blue dotted), 44 - 45 (green dashed), 45 — 46 
(yellow dot-dash), > 46 (red triple-dot-dash), compared to observations of Hasinger et al. (.2005.1 (filled circles of corresponding color, normalization adjusted 
for obscuration). 



peak luminosity. The n(Lpeak) distribution which produces the 
LDDE seen in Figure |3 has a simple functional form: a log- 
normal, with constant narrow width and normalization, and a 
center/peak directly related to the observed break L*(z) in the 
obse rved QLF. In keeping with the interpretation of the QLF 
from'Hopk ins et ai] j2005c). in which the faint-end is made up 
of brighter sources (dominated, of course, by the peak in the 
n(Lpeak) distribution at Lpeak ~ L^,(z)) in dimmer stages of their 
evolution, we have demonstrated that the observed LDDE is 
entirely accounted for by the quasar lifetime as a function of 
peak luminosity, not by any change in the shape of the peak 
luminosity or black hole mass distribution. 

A more detailed prediction of the faint-end slope should 
account for several points. Sample selection effects from 
reddening and extinction can be significant, especially if 
the fraction of obscured objects is a function of luminosity 
(e.g.,^eda et al. 2003; Hockins et al. 2005e). Our predictions 
agree best with observations from the hard X-ray, where the 
effect of attenuation is minimized, but a more accurate predic- 
tion of 7 in other wavebands must ac count for the joint evolu - 
tion of obscuration and luminosity jHopkins et alJl2005die^ ■ 
To probe very low accretion rates L <C 10~^ Lpeak, higher- 
resolution simulations, with more sophisticated models for 
low-efficiency accretion (rates well below Bondi or Edding- 
ton) and spectral modeling of the corresponding radiatively 
inefficient accretion flows are needed. Further, the contribu- 



tions to the faint-end QLF from quasars with Lpeak t^L* (see 
Hopkins et al. 2005e) will introduce overall curvature into the 
QLF (weakening the observed break, as observed by e.g. Wolf 
et al. 2003; Richards et al. 2005), which may aid in constrain- 
ing the «(Lpeak) distribution. Finally, at low enough L (typical 
of LlNERs and low-L Seyferts), we expect different, possi- 
bly stochastic, fueling mechanisms to contribute, which must 
be modeled to accurately predict the QLF shape below some 
limit. 

Still, our simple modeling of quasar feedback, and the novel 
picture of the QLF resulting from the application of real- 
istic, luminosity-dependent quasar lifetimes, accurately pre- 
dicts the faint-end slope of the QLF at all observed redshifts. 
Our modeling provides a simple, direct physical motivation 
for both the break luminosity L* and faint-end slope 7. The 
break L^, corresponds to the peak in the rate at which quasars 
with a given peak luminosity (final BH mass) are being cre- 
ated or activated at any given time, and 7 is determined (to 
first order) by the differential lifetime of these objects, as they 
spend substantial time at low L -C Lpeak, either in early stages 
of rapid BH growth to larger L or in sub-Eddington states in 
transition into or out of the brief period of peak quasar ac- 
tivity. The observed 7 and its evolution with redshift are a 
simple consequence of L* and its evolution. At high-z, L* 
is larger, implying most quasars have higher peak luminos- 
ity. From our simulations and analytical modeling of quasar 
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feedback, we expect higher-Lpeak objects to both grow and ex- 
pel gas more rapidly and violently, when they reach their final 
mass and feedback unbinds nearby material. Thus, brighter- 
Lpeak objects "die" more quickly, resulting in a flatter 7 as they 
spend less relative time in any given L < Lpeak- 

The observed values and evolution of 7 provide a direct 
test of the model of quasar lifetimes and the QLF proposed 
in Hopkins et al. (2005a-f), and are not predicted by mod- 
els which invoke simple "on/off" or pure exponential quasar 
light curves. Our model can predict 7 as a simple function 
of L,, itself tied to the characteristic quasars being created 
at any time, and is thus immediately useful for modeling the 
quasar contribution to reionization, which depends on 7 at 
high-z. For pu rposes where 7 is i mportant, the quasar life- 
times fitted in Hopkins et al. (2005e) should be slightly mod- 
ified to be Schechter functions with faint-end slopes 7(Lpeak) 
[the weak dependence of normalization on Lpeak found therein 
is entirely contained in 7(Lpeak)]- We have also provided sim- 
ple analytical forms for the quasar light curve, for use in semi- 
analytical models and other theoretical models of which re- 
quire the time-dependent quasar light curve (not simply the 
statistical properties contained in the quasar lifetime fits we 
have previously calculated), and which cannot resolve the de- 
tailed time dependence of quasar activity in individual merg- 
ers and interactions. We have demonstrated that the observed 
evolution of 7 and corresponding "luminosity-dependent den- 
sity evolution" of the quasar population are not just accounted 
for, but actually predicted in our model of the QLF, as a con- 
sequence of quasar feedback being more violent in higher 
peak luminosity (larger BH mass) systems, and the fact that 
the observed faint-end QLF is dominated by sources with in- 
trinsically brighter peak luminosity (Lpeak ~ L*(z)) in dimmer 
stages of their evolution. 

While more accurate predictions require additional detailed 
modeling, the results presented here allow future observations 
of 7 to directly constrain the differential quasar lifetime, even 
at L ^ Lpeak, as our modeling shows that the faint-end QLF 
effectively traces the quasar lifetime of Lpeak ^ L, sources, 
not the intrinsic source (e.g. peak L or BH mass) distribution. 
Such observations will further limit models of the distribution 
of quasar masses and host properties [through «(Lpeak)], and 
models of quasar fueling mechanisms and accretion (through 
the form of the lifetime/slope at low-L). As is clear in Fig- 
ure|3] the QLF faint-end slope as a function of redshift is still 
only poorly constrained by observations. Improved observa- 
tional constraints, in particular samples with a uniform selec- 
tion criteria which can span the range of redshifts z = 0-6, 
for which the selection criteria could be accurately modeled 
and compared in detail with observations, would provide very 
valuable constraints on theoretical models of quasar feedback, 
e.g. determining the relative contributions to the faint-end 
QLF from quasars growing to larger luminosities or relaxing 
after their peak, and constraining the efficiency of coupling 
between accretion energy and the surrounding ISM. Samples 
which cover the relevant redshift range in a uniform rest frame 
wavelength are also especially valuable, and measurements in 
different wavelengths can provide different constraints. For 
example, the faint-end hard X-ray QLF may be more robust 
against obscuration effects and so provide a better indicator of 
the bolometric luminosity function, tracing e.g. relatively low- 
luminosity stages hidden by circumnuclear starbusts. The op- 
tical QLF, on the other hand, is in our modeling more closely 
associated with the peak quasar luminosity and "blowout" 
phase, meaning that the faint-end optical QLF (measurement 



of which would require extending the completeness of current 
high-redshift optical surveys by several magnitudes) can pro- 
vide an valuable constraint on the underlying peak luminosity 
distribution, a particularly valuable quantity in itself because 
it determines, in our modeling, the QLF in all other bands, and 
reflects (as well as constrains) much more directly the under- 
lying cosmological context, such as merger rates as a function 
of host galaxy mass. The combination of the two, yielding 
a reliable bolometric faint-end luminosity function slope and 
underlying peak luminosity distribution, would allow mea- 
surements of 7(z) to be translated reliably into 7(Lpeak), and 
could be de-convolved to construct an entirely observational 
determination of the quasar lifetime as a function of luminos- 
ity. These observations, though difficult at very high z, can 
probe redshifts at which direct observations of quasar hosts 
and masses are inaccessible, constraining theories of early 
quasar evolution. 

Combined with the modeling presented in Hopkins et al. 
(2005a-f), this motivates a completely self-consistent picture 
of the quasar lifetime and luminosity function, derived from 
the input physics of our simulations and without arbitrary tun- 
able parameters. With the critical recognition that the quasar 
lifetime is luminosity-dependent, with quasars spending more 
time at low luminosities than their peak luminosities (as has 
now also been suggested by observations, e.g. Adelberger & 
Steidel 2005), we have shown that a large range of observed 
properties of the quasar and spheroid population are predicted 
to high accuracy, over a wide range of observed redshifts. 
These works also provide several simple tests of our model 
of the quasar lifetime, and correspondingly the model of the 
faint-end slope presented here. In Hopkins et al. (2005a-e), 
we explicitly demonstrate that these predictions include, e.g. 
the QLF at frequencies from optical to hard X-ray and at many 
observed redshifts, the distribution of Eddington ratios as a 
function of luminosity and redshift, the column density dis- 
tribution at various wavelengths, the fraction of broad-line 
quasars as a function of luminosity and redshift, the active 
BH mass function of both Type 1 and Type II quasars, the 
relic BH mass function and total density, the cosmic X-ray 
background spectrum, the anti-hierarchical evolution of BH 
mass, the radio source population at high redshift, correla- 
tions between IR emission, star formation, and quasar ob- 
scuration, and the quasar lifetime as a function of luminosity 
and host galaxy properties. The distribution of Eddington ra- 
tios as a function of luminosity is a direct test of this model, 
and the degree to which typical Eddington ratios decrease be- 
low the break in the observed luminosity function informs the 
extent to which the faint-end slope is dominated by sources 
at Lpeak ^ L,. Likewise, the active BH mass function con- 
strains the n(Lpeak) distribution, and an accurate measurement 
of its peak can be used to determine the Lpeak which dominates 
the faint-end slope, providing a direct means to constrain the 
slope 7 as a function of Lpeak- The broad line fraction as a 
function of luminosity can provide a similar constraint, as the 
luminosity at which the QLF transitions to being broad-line 
dominated is of order the turnover luminosity in the Schechter 
function fit to the quasar lifetime for Lpeak ^ the peak of the 
n(Lpeak) distribution. 

Further, Lidz et al. (2005) show this model also predicts the 
correlation function and bias of quasars as a function of lumi- 
nosity and redshift, and the quasar-galaxy cross-correlation 
function, dramatically different from the prediction in models 
with a simple quasar lifetime. The mean quasar bias (not- 
ing that the bias depends on luminosity much more weakly in 
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our model than in simpler models of the quasar lifetime) is a 
direct probe of typical quasar host halo masses, constraining 
the peak in the n(Lpeak) distribution and thus allowing for a 
determination of 7(Lpeak)- The steepness of the bias vs. lu- 
minosity is sensitive to the width of the n(Lpeak) distribution, 
and thus limits the contribution of quasars with Lpeak "4^* to 
the faint-end slope. Thus far the observations (e.g. Croom et 
al. 2005, Adelberger & Steidel 2005) support a picture of a 
peaked n(Lpeak) distribution in close agreement with our mod- 
eling, but the bias as a function of luminosity is still only 
poorly determined by observations. 

Hopkins et al. (2005f) extend these predictions to the prop- 
erties of spheroids and red galaxies, self-consistently and 
accurately predicting (from the observed quasar luminosity 
function) the distribution of velocity dispersions, the mass 
function, mass density, and star formation rates, the luminos- 
ity function in a wide range of observed wavebands, the total 
red galaxy number density and luminosity density, the dis- 
tribution of colors as a function of magnitude and velocity 
dispersion for various color-magnitude combinations, and the 
distribution of formation redshifts/ages as a function of both 
mass and luminosity, each as a function of redshift at different 
observed redshifts. Our modeling also reproduces the X-ray 
emitting gas properties (Cox et al. 2005), morphologies and 
metallicities (Cox et al. 2005, in preparation), Mbh - o" rela- 
tion (Di Matteo et al. 2005) and fundamental plane relation 
(Robertson et al. 2005, in preparation), the dispersion about 
each relation (Robertson et al. 2005), and bimodal color dis- 
tribution (Springel et al. 2005a) of observed galaxies. The 
distribution of these properties and e.g. spheroid luminosity 
and mass functions can be used to constrain the form of the 
quasar lifetime based on purely observational arguments, by 
using observed correlations between e.g. stellar or virial mass 
and BH mass in the context of the merger hypothesis of ellipti- 
cal galaxy formation to determine the underlying rate at which 



galaxies of a given final (combined) stellar mass must merge 
to produce the remnant spheroid population and its evolution 
with redshift (essentially determining the differential growth 
with time and mass) and from this the corresponding n(Lpeak) 
distribution. Having determined the rate at which quasars 
of a given peak luminosity are created, a de-convolution of 
the observed quasar luminosity function and «(Lpeak) can con- 
strain the allowed range of quasar lifetimes as a function of 
peak luminosity (essentially the reverse of our modeling, in 
which we de-convolve the QLF and quasar lifetime to deter- 
mine n(Lpeak))- 

The wide range of observations predicted by this model 
provides a number of very different means by which to ob- 
servationally test both the qualitative picture proposed and 
the detailed dependence of quasar lifetime and correspond- 
ing faint-end QLF slope on quasar peak luminosity. We have 
demonstrated that a proper modeling of quasar evolution, in- 
cluding the critical effects of quasar feedback on the host 
galaxy in major mergers, can completely account for (and 
in fact predicts) the detailed evolution of the faint-end QLF 
slope and luminosity-dependent density evolution, without fit- 
ting to the faint-end observations or fine-tuning any parame- 
ters, and is entirely self-consistent with the large number of 
detailed observations above, as a simple consequence of real- 
istic, physically motivated, luminosity-dependent quasar life- 
times. 
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